2. 50 observaciones de una \(\mathcal{N}(0,1)\) y otras 50 observaciones de una \(\mathcal{N}(2,1)\). ¿Cómo se verán las 100 caras de Chernoff si \(_X_\) y \(_Y_\) definen la lÃnea de la cara y la oscuridad del cabello?, ¿esperan caras similares?, ¿cuántas caras lucen como oservaciones de \(_Y_\) cuando aún son observaciones de \(_X_\)?
ind = matrix(0, ncol = 36) # define una indicadora para el argumento which.row
ind[,13] = 1 # linea derecha de la cara
ind[,14] = 1 # oscuridad del cabello lado derecho
ind[,31] = 1 # linea izquierda de la cara
ind[,32] = 1 # oscuridad izquierda del cabello
x = rnorm(50)
y = rnorm(50, mean = 2)
z = t(cbind(t(x),t(y))); # arma matriz (100x1)
faces(as.matrix(z[1:50]),ind, main="Observaciones 1 a 50", ncol.plot = 5, nrow.plot = 10) # primeras 50 caras
## effect of variables:
## modified item Var
## "height of face " "Var1"
## "width of face " "Var1"
## "structure of face" "Var1"
## "height of mouth " "Var1"
## "width of mouth " "Var1"
## "smiling " "Var1"
## "height of eyes " "Var1"
## "width of eyes " "Var1"
## "height of hair " "Var1"
## "width of hair " "Var1"
## "style of hair " "Var1"
## "height of nose " "Var1"
## "width of nose " "Var1"
## "width of ear " "Var1"
## "height of ear " "Var1"
faces(as.matrix(z[51:100]),ind, main="Observaciones 51 a 100", ncol.plot = 5, nrow.plot = 10) # primeras 50 caras
## effect of variables:
## modified item Var
## "height of face " "Var1"
## "width of face " "Var1"
## "structure of face" "Var1"
## "height of mouth " "Var1"
## "width of mouth " "Var1"
## "smiling " "Var1"
## "height of eyes " "Var1"
## "width of eyes " "Var1"
## "height of hair " "Var1"
## "width of hair " "Var1"
## "style of hair " "Var1"
## "height of nose " "Var1"
## "width of nose " "Var1"
## "width of ear " "Var1"
## "height of ear " "Var1"
x
## [1] -0.63709002 -1.55517217 0.97120839 0.25899932 0.62262458
## [6] 0.88735305 0.84519325 -1.45545610 0.18242179 -0.38505343
## [11] -0.29708418 -0.87640833 1.94051419 1.61181490 1.24085868
## [16] -1.14741563 -0.13287059 0.15123131 -1.65420017 -0.61784681
## [21] 0.61044777 -0.63901666 0.44668333 1.60617987 0.80826290
## [26] 0.54959329 0.78147470 0.40767688 -1.92743791 -1.63081938
## [31] -1.46981445 -0.43132915 0.03896066 -1.35444287 0.22807248
## [36] -0.88071871 0.33247103 0.77120602 -0.48212797 -1.14623617
## [41] 1.02019694 0.42310914 0.80433054 0.76080547 0.56774406
## [46] 0.85727204 -0.21988876 -0.79397466 0.18934862 -1.05086577
y
## [1] 3.0882151 2.2798873 2.2752659 4.0571267 0.8423817 2.6949473
## [7] 2.3832699 1.2253625 -0.2155235 0.8185916 2.4536493 2.0031057
## [13] 2.9490554 2.5933593 2.2162759 3.4064165 1.9195255 3.0370970
## [19] 2.2699174 1.5282019 -0.6730547 2.5954333 3.0287296 2.2045015
## [25] 0.7375079 1.7323372 2.2546905 0.5836871 2.0513613 4.2506394
## [31] 1.5296861 0.4044310 2.0820663 1.8680234 2.5995023 2.5387882
## [37] 2.6384141 4.3536440 2.6024578 2.6840687 2.6300256 2.5174938
## [43] 2.6944678 2.9935385 3.1770656 1.3269875 0.6731204 3.4430816
## [49] 2.6288425 1.9472525
Sà hubo caras similares. Hay 4 payasos en Y comparado con los 3 que hay de X. Hay casi la misma cantidad de caras rojas y cabello verde; hay 4 caras amarillas en Y comparadas con las 10 que hay en X. Las caras rojas son aquellos valores cercanos a 0, mientras que los payasos son los valores más grandes (+/- 3)
3. Consideren los siguientes datos
a. Encuentra la proyección de \(X_1\) sobre 1’=(1,1,1,1,1,1)
# x1 = Nómina de jugadores
x1 <- c(3497900,2485475,1782875,1725450,1645575,1469800)
# Defino al vector de 1 normalizado
ones <- rep(1,6)/sqrt(6)
# Definición de proyección de x1 sobre 1
proy = t(x1) %*% ones %*% ones
proy
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 2101179 2101179 2101179 2101179 2101179 2101179
b. Calcula el vector desviación. Relaciona su longitud a la desviación estándar.
# Definición de vector desviación
vDesv <- x1 - proy
vDesv
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 1396721 384295.8 -318304.2 -375729.2 -455604.2 -631379.2
# Desviación estándar
StanDev <- sd(x1)
StanDev
## [1] 767752.2
Notemos que el cuadrado de la longitud del vector desviación es igual a \(n\times d^2\) donde d es la desviación estándar
# Longitud del vector
norma <- norm(vDesv)
norma^2
## [1] 1.950829e+12
6*(StanDev^2)
## [1] 3.536661e+12
c. Grafica (a escala) el triángulo formado por \(y_1\), \(\bar{x}_1\), \(y_1-\bar{x}_1\). Identifica la longitud de cada vector en tu gráfica
nx1 = norm(x1, type = "2")
nproy = norm(proy, type = "2")
nvDesv = norm(vDesv, type = "2")
x = c(0, nproy/nx1, nproy/nx1)
y = c(0, 0,nvDesv/nx1)
plot(x, y, xlim = c(0,1.1), ylim = c(0, .4))
arrows(0,0, x1 = nproy/nx1, y1 = 0)
arrows(0,0, x1 = nproy/nx1, y1 = nvDesv/nx1)
arrows(nproy/nx1,0, x1 = nproy/nx1, y1 = nvDesv/nx1)
vectores = cbind(x,y)
text(vectores, labels = c(round(nproy/nx1,3),round(nvDesv/nx1,3),1), adj = c(0,-1))
d. Repetir los incisos (a) a (c) para \(X_2\) a’. Encuentra la proyección de \(X_2\) sobre 1’=(1,1,1,1,1,1)
# x2 = % de perdidos-ganados
x2 <- c(0.623,0.593,0.512,0.5,0.463,0.395)
# Definición de proyección de x2 sobre 1
proy2 = t(x2) %*% ones %*% ones
proy2
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.5143333 0.5143333 0.5143333 0.5143333 0.5143333 0.5143333
b’. Calcula el vector de desviación
# Definición de vector desviación
vDesv2 <- x2 - proy2
vDesv2
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.1086667 0.07866667 -0.002333333 -0.01433333 -0.05133333 -0.1193333
# Desviación estándar
StanDev2 <- sd(x2)
StanDev2
## [1] 0.08376555
c’. Grafica (a escala) el triángulo formado por \(y_2\), \(\bar{x}_2\), \(y_2-\bar{x}_2\). Identifica la longitud de cada vector en tu gráfica
nx2 = norm(x2, type = "2")
nproy2 = norm(proy2, type = "2")
nvDesv2 = norm(vDesv2, type = "2")
x = c(0, nproy2/nx2, nproy2/nx2)
y = c(0, 0,nvDesv2/nx2)
plot(x, y, xlim = c(0,1.1), ylim = c(0, .4))
arrows(0,0, x1 = nproy2/nx2, y1 = 0)
arrows(0,0, x1 = nproy2/nx2, y1 = nvDesv2/nx2)
arrows(nproy2/nx2,0, x1 = nproy2/nx2, y1 = nvDesv2/nx2)
vectores = cbind(x,y)
text(vectores, labels = c(round(nproy2/nx2,3),round(nvDesv2/nx2,3),1), adj = c(0,-1))
e. Grafica (a escala) los dos vectores desviación \(y_1-\bar{x}_1\) y \(y_2-\bar{x}_2\). Calcula el valor del ángulo entre ellos
f. Calcula la varianza muestral generalizada (vmg) det(S) para estos datos e interpreta
X <- cbind(x1,x2) # Creo la matriz de X1 y X2
n <- nrow(X)
uno <- rep(1,n)
xBar <- t(X) %*% uno / n # x barra
H <- diag(n) - uno%*%t(uno) / n # H
Sn <- t(X) %*% H %*% X / n # Sn
S <- n/(n-1) * Sn
vmg <- det(S)
vmg
## [1] 844182191
Como la vmg es proporcional al cuadrado del olumen generado por los p vectores de desviación. COmo vmg = 844182191, entonces sabemos que el volumen del elipsoide al cuadrado generado por \(S\), será igual al vmg multiplicado por una constante proveniente de la ecuación del elipsoide correspondiente a \(S\).
g. Calcula la varianza muestral total (vmt) tr(S) para estos datos e interpreta
vmt <- sum(diag(S))
vmt
## [1] 589443426354
Geométricamente, la vmt es la suma de las longitudes al cuadrado de los p vectores de desviación divididos entre n-1; sin embargo, no le presta atención a la orientación de los vectores residuales. Como vmt = 589443426354, entonces sabemos que estas corresponde a la suma de las varrianzas de los p elementos de la matriz.
4. Dibuja las elipsoides sólidas para las tres matrices siguientes y determina los valores de los ejes mayores y menores
5. Archivo del INEGI
a. Configura la matriz X para poder operar con ella
# datos <- read.csv("D:\\ITAM\\Aplicada III\\Tareas\\HW2\\INEGIConstruccion2017.csv", header=FALSE)
datos <- read.csv(".\\INEGIConstruccion2017.csv", header=FALSE)
datos <- datos[-1,]
names(datos) <- as.matrix(datos[1,]) # Tomo primer renglón
datos <- datos[-1,] # Elimino el primer renglón
datos[] <- lapply(datos, function(x) type.convert(as.character(x))) # Primera columna names(datos) se vuelve el header
# Datos de enero 2017
datos <- datos %>% filter(Periodo == "2017/01")
head(datos)
datosTotalesNac <- datos[,str_detect(names(datos), "Total nacional")]
# Horas trabajadas (dependiente + no dependiente, o lo mismo que es la suma de los totales por estado)
hrasTrabajadas <- datosTotalesNac %>% select(2,6) #dependiente y no dependiente
# x1: total de horas trabajadas
sumaHorasT <- hrasTrabajadas[,1]+hrasTrabajadas[,2]
#x2: valor total de producción
totalProd <- datos[,str_detect(names(datos), "tipo de obra > Total Total")]
#x3: total de horas trabajadas Dependiente
horasDep <- hrasTrabajadas[,1]
#x4: Obreros Dependiente
obrerosDep <- datosTotalesNac[,str_detect(names(datosTotalesNac), "Obreros")]
#x5: Empleados Dependiente
empDep <- datosTotalesNac[,str_detect(names(datosTotalesNac), "Empleados")]
#x6: Propietarios Dependiente
propDep <- datosTotalesNac[,str_detect(names(datosTotalesNac), "Propietarios")]
#x7: total de horas trabajadas No dependiente
horasNoDep <- hrasTrabajadas[,2]
(X <- c(sumaHorasT, totalProd, horasDep, obrerosDep, empDep, propDep, horasNoDep))
## [1] 129460.750 33007656.590 107742.613 84562.852 21536.763
## [6] 1642.998 21718.137
b. Calcula el vector de medias y la matriz de covarianzas de la matrix X
#Si consideramos la X como descrita arriba, el vector de medias es el mismo que el vector X